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Abstract 



We study the canonical solution of a family of classical n — vector spin 
models on a generic (i-dimensional lattice; the couplings between two spins 
decay as the inverse of their distance raised to the power a, with a < 
d. The control of the thermodynamic limit requires the introduction of a 
rescaling factor in the potential energy, which makes the model extensive 
but not additive. A detailed analysis of the asymptotic spectral properties 
of the matrix of couplings was necessary to justify the saddle point method 
applied to the integration of functions depending on a diverging number 
of variables. The properties of a class of functions related to the modified 
Bessel functions had to be investigated. For given n, and for any a, d and 
lattice geometry, the solution is equivalent to that of the a = model, 
where the dimensionality d and the geometry of the lattice are irrelevant. 



Key words: Classical n — vector models; long-range interactions; mean 
field solution; ensemble equivalence. 



1 Introduction 



Systems with long-range interactions do not have a well defined thermody- 
namic limit/ 1 ) and their equilibrium statistical mechanics is not inside the 
framework of the more studied short-range systems. For a pair interaction 
long-range means, in this framework, that the modulus of the potential 
energy decays, at large distance, not faster than the inverse of the dis- 
tance raised to the power of the spatial dimension. In fact, conditions 
can be given for the existence of the thermodynamic limit/ 2 ' 3 - 1 one being 
temperedness, that for pair interactions is exactly the requirement that the 
coupling is not long-range as just defined, and the other being stability, 
which requires that the global minimum of the potential energy does not 
diverge to — oo more than linearly with the number of degrees of freedom. ( 4 ) 
The physically relevant example of continuous long-range systems violating 
temperedness is that of particles interacting with gravitational forces; here 
there is in addition the problem of the potential energy at short distance, 
that manifestly violates stability. Both conditions are necessary for the 
system to be extensive and additive. Extensivity means that the specific 
thermodynamic potentials (i.e., the thermodynamic potentials per parti- 
cle), do not diverge in the thermodynamic limit, while additivity means 
that the thermdynamic potentials of the whole system are, in that limit, 
the sum of those of its component macroscopic parts. ( 5 ' 6 - ) In particular, 
additivity is essential if one wants to derive the canonical ensemble from 
the microcanonical ensemble. 

Since the classic paper by Ising/ 7 -* magnetic models on a lattice are the 
most used to investigate the statistical physics of interacting many-body 
systems, because their mathematical tratment, although still difficult, is 
often more affordable than that of continuous systems. Many generaliza- 
tions of the Ising model have been considered, over time; among them, also 
models with long-range interactions. In lattice systems the problem of the 
behaviour of the potential energy at short distance is not present, and no 
collapse, as in gravitational models, can exist; the equality of all the co- 
ordinates, as in the case of maximum magnetization, is not a pathological 
configuration. The lack of stability and temperedness can be related only 
to a coupling between spins decaying too slowly with distance. Thus, the 
study of long-range magnetic systems on lattices can give insights on the 
statistical properties in general of long-range potentials, and is therefore 
useful to investigate how to extend the framework of statistical mechanics 
to systems that do not obey stability and temperedness. 

In this paper we study the equilibrium statistical mechanics of classical 
n — vector spins (i.e., n-dimensional unit vectors) fixed on the nodes of 
a generic rf-dimcnsional lattice, and interacting pairwise through a long- 
range potential. We consider generic values of n and d. The couplings that 
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we consider decay as the inverse of the distances between spins raised to 
the power a; if a is not larger than d the interaction is not tempered, and 
the specific energy (or energy density) of the system diverges in the ther- 
modynamic limit. This divergence can be cured with a Kac's prescription, 
that is gauging the potential energy with an appropriate scaling function 
of N, the number of spins, and dS 8 ^ Then extensivity is enforced through 
a control of the thermodynamic limit but, due to the long-range couplings, 
additivity does not hold, and ensemble equivalence, whose proof is based 
on the possibility of separating the energy of a subsystem from that of the 
whole, might not be guaranteed in that limit/ 6 ) We study our system in 
the canonical ensemble and we find the exact soultion, and in subsection 
4.5 we will show that for this class of systems the microcanonical and the 
canonical ensembles are equivalent, in spite of the non additivity. 

The main motivation of this work is to show that all the systems in- 
cluded in our study, and obtained for the different values of n, d and for 
the different possible geometries of the lattice, are essentially all equiva- 
lent as far as their equilibrium statistical thermodynamics is concerned; 
some quantitative differences (e.g., the value of the critical temperature) 
are found between systems with different n, but the overall behaviour is 
always the same. This is true also with regard to the dependence on the 
distance of the coupling between spins, as long as it remains long-range: 
although we consider a power law function, we argue in the final discussion 
that any other form of the long-range coupling would have brought us to 
the same results. 

The paper is divided into five sections. In Section 2 we present the class 
of long-range magnetic models that we consider. The canonical solution 
for the mean field case (a = 0) is given in Section 3. Section 4 is devoted to 
the general mean field case (0 < a < d); it contains the main results of this 
paper and is divided into five subsections: in the first three the canonical 
solution is obtained, while the last two are devoted to the justification of 
the saddle point method and to the proof of the ensemble equivalence. In 
the discussion of Section 5 we argue about possible generalizations of our 
results. 



2 The Model 

We consider here classical models which belong to the general class of 
the n — vector models; ^ in our previous work^ 10 - 1 we had studied a kind 
of XY-model, which is a vector model with n = 2; the case n = 3 would 
correspond to a classical Heisenberg model. Classical n— vector spin models 
correspond to the "infinite spin" limit of magnetic quantum systems;^ 11 -' 
these models are related also to the spherical model of Berlin and Kac/ 12 ) 



2 



in the limit n — > oo the n — vector model reduces to the spherical model/ 13 - 1 
The consideration of long-range couplings in magnetic models on a lattice 
was started in the sixties of the last century, when the basic mathematical 
techniques were established. *- 14 ) 
The hamiltonian we consider is: 



N 

I 

H 



1 N 

- £ J«(l - Si ■ Sj). (1) 

For each i the spin Si is a unit n-dimensional vector; its position can 
be specified by n — 1 angles and its cartesian coordinates are related to 
these angles through the definition of the polar coordinates in lZ n . The 
index i = 1, . . . , N labels the sites of a generic d-dimensional lattice, with 
d integer. We conventionally put an external factor 1/2 and assume the 
rescaled form 1 — Si ■ Sj which allows a free choice of the diagonal terms 
J a because of the constraint Sf = 1. We consider interactions in the form 
of an inverse power of the distance between lattice site i and j: 

Jij = ^ (2) 

with periodic boundary conditions and the nearest image convention for the 
distance r.y; a > sets the range of the interaction, which is long-range if 
a < d and short-range if a > dM 5 ' Different models remain characterized 
by different values of d, n and a. As said above, most of the work has been 
done in the study of models with a > d; a = oo is the limit of nearest 
neighbors interactions. We focus our attention instead on the case a < d\ 
in this case the interaction (Q) causes the energy (|l|) to be not extensive; 
thus the partition function does not admit a well defined thermodynamic 
limit. The problem can be solved through the rescaling: 



N, 



■I, ■ ^ (3) 



where N is a function of N, a, d and the geometry of the lattice. The 
above substitution is a generalization of the usual one found in the mean 
field version of this kind of models: the general interaction Jy is replaced 
by a constant one (i.e., a = 0), but a rescaling factor 1/N is necessary to 
control the thermodynamic limit. The mean field substitution is recovered 
as a particular case of (||) for a = 0. In the following we will sometimes 
refer to the general model, with < a < d, as the general mean field case. 

We previously studied this rescaled model in the case n — 2;( 10 ) in 
that short paper, where only few mathematical arguments and details were 
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given, we therefore considered only the model represented by planar rota- 
tors on a lattice, but with the addition of a kinetic term in canonical vari- 
ables, conjugate to the angles. Here we will not consider the kinetic energy, 
that in the canonical ensemble gives an additive trivial contribution to the 
thermodynamic potentials; we will make very short comments on this in 
the right places. Computation of the partition function for n = 2 has shown 
its universality in a: the free energy does not depend on the value of a < d 
and is thus equal to the mean field one for a — 0. The model of planar 
rotators, n = 2, has also been the subject of some studies on the dynamical 
properties: the universality of the thermodynamics in a one-dimensional 
lattice (d = 1) with respect to a < 1 was suggested by the numerical study 
in ref. 16; interesting metastable states have been observed for a — O^ 17 ' 18 - 1 
and for a < 1 with d = lS ls ' 

Here we take a general n — vector model and we show how to treat 
the interaction 1 /rg in the case a < d as far as the computation of the 
partition function is concerned. 



3 The mean field case a = 



In this Section we review the canonical solution of the a = model; in the 
case of planar rotators, that is n = 2, the solution has already been given in 
ref. 20. For a = the interaction is homogeneous, i.e., each spin interacts 
with all the others with the same strength. In this case the rescaling (||) is 
usually performed with N = N. The hamiltonian is then: 



N N 

H = w Ea-Si-s.o-h-^Si, (4) 

where we also include an external magnetic field h, which is coupled in the 
usual way to the spin vectors Sj. The partition function for the model is 

Z = Jd N Oe-> 3H , (5) 

where (3 is the inverse temperature, dQi is the surface element of the unit 
sphere in dimension n > 2, and d N <d = dQi . . .dO^; and for n = 1 the 
integral is replaced by a sum on all the possible Ising spin configurations: 
Esi=±i s N =±i e _/3ff - Here, and later for the general case, we do not 
consider the kinetic part, that in the classical partition function trivially 
decouples. The magnetization is given by: 

Me-(Vs,) = Tr-InZ, (6) 

N x <^ ' Nf3dh ' w 
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where (•) = J d N Q(-)e 1311 is the usual canonical average. With B = /3h 
and C = exp(— N/3/2), we can rewrite: 

Z = C Jd N Oe^M 2 +*-E> s '. ( 7 ) 
Using the gaussian transformation 

1 r+°° 2 

e aS ~ = -jL= / dze~^ +Sz a > (8) 



on each term of the square modulus of the vector Si (we emphasize that 
the above expression is valid when a > 0, or more generally, for complex 
a, when its real part is greater than 0), we linearize the quadratic term in 
(0) and obtain: 

Z = C (JL} ^ I dz e-^ 2 J d N e £< S '-( z+B >, (9) 

where dz = dzi . . . dz n . Here and in the following the notation b 2 will 
denote, interchangeably with b 2 , the scalar product of the vector b with 
itself, i.e., its square modulus. The last integral separates on the sites i 
and gives N identical contributions, the functional form of which depends 
on the spin dimension n. In Appendix A we show a more convenient way 
to write these integrals on the unit sphere and how they are related to 
the modified Bessel functions; besides, we introduce the notation, that 
expresses the surface integrals in ([)]) in terms (for proper values of n) of 
a function G n (x) and of the area f2„ of the unit sphere in n dimensions. 
Some properties of the functions defined in Appendix A are needed for 
the analysis of this mean field case and of the more general case with 
a ^ 0; these properties are proved in Appendix B. Following the notation 
introduced above we rewrite the partition function as: 

i ( JL\ 1 ' 2 [ <fafi ^[-^- a +taG„_ a (|-+B|)] 



Z = CK-i [2^) J dz e [ ~ W * +mw - 2UZ+0|j J. (10) 

The integral is computed with the saddle point method; therefore, we need 
to find the stationary points of the function in square brackets in the ex- 
ponent, and consider those that are maxima; the dominant contribution to 
the integral will be determined by the absolute maximum. If we call /(z) 
the function in square brackets in (|l^), the stationary points are given by 
the solutions of the system of n equations (one for each component of 
z): 

-f =- ^ +9n _ 2 (|z + B|)f^-f =0, M = l,2,...,n (11) 
oz^ [j l z + B| 
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where the function g is the logarithmic derivative of G. The free energy 
per particle (or specific free energy) will be given by: 



8F = lim — bxZ 

N^oo N 

= -- + lnf2 n _ : 



_J_ z 2 +lnGn _ 2 (| z + B | 

2p 



= -^ + lnn n _i +max/(z). (12) 

We note that the hessian of /(z) in the maximum does not appear in (|l2|), 
since its contribution becomes vanishingly small in the thermodymamic 
limit, N — > oo. An analogous argument can be used for possible degenera- 
cies of the absolute maximum (see few lines below). The study of ([ll]) is 
presented in some details in Appendix C; here we only show the results and 
the final expressions for the magnetization and for the equation of state. 
Eq. ( pi] ) can have more than one solution, depending on the value of 8 
and B = |B|. In any case, the relevant stationary point z* is such that its 
modulus z* satisfies the self-consistency equation: 

z* = Pg„- 2 (z* + B), (13) 

which is a generalization of the Curie- Weiss equation found in the mean- 
field solution of the Ising model. In fact, as shown in Appendix B, g n (x) 
has the same qualitative features of tanh(x). A visual aid (Fig. 1) is 
provided in Appendix C together with the details of the study of Eq. ( |TT| ) . 
To complete the solution, we must add the following specifications. When 
B > (i.e., h = |h| > 0), for which Jl3| ) has a positive solution z* > 0, z* 
is parallel to h. When h = 0, we have to distinguish between 8 > B c = n 
and 8 < Be'- for inverse temperatures not greater than the critical value 
8 C the only solution of (|l3| ) for h = Q is z* = 0; instead, for 8 > 8 C 
there is also a positive solution, and this is the relevant one. In this last 
case, the direction of z* (if n > 1) remains undetermined; in other words, 
the stationary point is infinitely degenerate, or doubly degenerate if n = 1. 
Therefore, to be more precise, one should then perform in ([l^) , when 8 > 8 C 
and h = 0, an integration over the angular coordinates of z (or a sum over 
the two directions if n = 1) before applying the saddle point. This would 
give in ( |l2[ ) a further factor ln 3> , which does not contribute to F in the 
thermodynamic limit. 



The magnetization M is given by (|6|), that applied to (12) gives M = 
Since in the case h = and 8 > 8 C the direction of z* is not determined, 
because of degeneracy, the actual direction of M in a real system is deter- 
mined by the boundary conditions, and there is a spontaneous simmetry 
breaking. The magnetization modulus M = |M| becomes zero at f3 = 8 C 
and remains zero for 8 < 8 C ; at 8 C there is a second order phase transition. 
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Finally, the equation of state, relating the specific (potential) energy 
U = jf(H) to the temperature (through the magnetization modulus M), 
is given by: 

U = - lim J-A] n ^=I(i_ m 2 ) - hM. (14) 
n^oo Nd/3 2 v ; y ' 

If we had considered also the kinetic energy, then in (JlJ) we would have 

had an additional term ^y^mf^^, an( l consequently in ( |l4|) a further 

term would have appeared, making in that case U the total specific 
energy. For (3 < j3 c the specific potential energy remains constant equal 
to i (since M remains equal to zero), and only the specific kinetic energy 
increases. 

The phase space volume at disposition of the rotators increases with n, 
and this is reflected in a critical temperature T c — — decreasing with n; 
apart from this quantitative difference, the overall behaviour is the same 
for all values of n. 



4 The general mean field case < a < d 

This is the central Section of the paper. We take into consideration the 
general form of the interaction, i.e. Jy — l/rfj with < a < d (see (Q) 
and (|J)). From the form of H it is clear, since Si is a unit vector, that the 
values of Ju can be chosen arbitrarily, as long as they are finite. We will 
use this freedom below, for the computation of the partition function. 

We do not explicitly consider the marginal case a — d. However, in 
the discussion we will comment on this point and on possible more general 
forms of the long-range couplings. 



4.1 Thermodynamic limit 

We have to consider the problem of the rescaling of the interaction param- 
eters, to control the thermodynamic limit. For the classical lattice systems 
considered here the existence of the thermodynamic limit is guaranteed by 
a restriction on the long-range part of the potential energy. For transla- 
tionally invariant interactions, i.e., when in our case Jy depends only on 
the distance vector from site i to site j, the restriction takes the form: 

N 

lim \ J ij\ < °°: ( 15 ) 

where i can take any value, since the translational invariance implies that 
the above sum is independent of i. It corresponds to an extensivity require- 
ment for the energy H. The rigorous demonstration of the sufficiency of 
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this condition for the existence of the thermodynamic limit of the partition 
function can be found in ref. 2, where the Ising model (n=l) is considered, 
but the procedure can be extended to the general case. 

For our system we introduce translational invariance also for finite N, 
through the use of periodic boundary conditions and a choice of the same 
finite value b for all the diagonal terms Ja\ this is convenient for many 
steps of the analysis of our model. Apart from this, condition ([l5|) requires 
the analysis of the quantity 

* = E^ (16) 

in the N — > 00 limit. We have dropped the j = i term, since the single 
value of J a does not determine the convergence property of (|l6|). To sec 
the large N behavior of ( |l6| ) for a < d we can shift to the integral, that 
gives S ~ N x - a / d (see, e.g., next subsection). Divergence occurs and 
must be compensated to have a thermodynamically well defined model. 
Compensation is readily obtained with the position (^), taking 

N ~ S. (17) 

4.2 Interaction matrix (spectral properties) 

The nonrescaled couplings Jij are the entries of a real, symmetric N x N 
matrix. In this subsection we study its eigenvalues, since this is important 
for the computation of the partition function. This matrix can be diag- 
onalized through a unitary matrix V and its eigenvalues are all real. In 
our computation of the partition function it will be important that all the 
eigenvalues are positive, since we will make use of the gaussian transfor- 
mation (^). We show here how we can use the freedom on the value of J a 
for this purpose. 

Let us denote the position of lattice site i by r^. We stress that in this 
subsection, and only in this one, boldface characters denote d-dimensional 
vectors of the lattice space or of its dual; in the rest of the paper, before 
and after this subsection, they denote n-dimensional vectors related to 
the dimensionality of the spins. The function is the distance between 
the point r» and the nearest image of r^, and then it is invariant under 
translations. The same is true for Jy, once we use the freedom on the 
values of Ja taking all these diagonal elements equal to the same constant 
b. Therefore, if we let = r j — r \ , where is the image of which is 
nearest to , denoting?^- = |ry| and introducing the notation Jy = J(rij), 
we have: 

J{r ij ) = { -i if ^ = ° (18) 
I r i] otherwise. v ' 
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Following these definitions, the eigenvalues of the matrix J are obtained 
through the d-dimensional Fourier transform 



A k = £ J(r)e~ ik r = J2 r- a e~* \ (19) 

r r 

where the sum is on all the lattice points and k denotes any of the N 
reciprocal lattice vectors contained in the first Brillouin zone. The reality 
of the eigenvalues follows from (|l9| ) and the properties of J(rij). It is also 
evident that Ao is the largest eigenvalue. If we isolate the r = term we 
have: 

A k = fe + ^r- Q e- lkr , (20) 

r^O 

that shows that the whole spectrum is linearly translated by b. The re- 
maining sum for k = corresponds to the sum S defined in (|l6|). It is clear 
that for a > d all the eigenvalues are finite in the thermodynamic limit. 
We now restrict to the case a < d. The large N behaviour of S can be 
estimated shifting S to an integral 

i 

pN~d 



If we take 



S = 2jr-"~ / r d - l - a dr ~ JV 1 ""*. (21) 



JV 

N = \ =b + S = J2 J ij, ( 22 ) 



the requirement ([17]) is satisfied as an equality in the thermodynamic limit 
because b is a finite quantity; as noted before, the last expression in ( p2[ ) 
does not depend on i. It is also possible to estimate the behavior of Ak 
for k ^ for large N. It is easy to see that the sum in (|20|), if k = |k| 
is different from 0, remains finite in the thermodynamic limit, and the 
behavior in k, again shifting to an integral, can be found to be 

This expression does not consider the sign of the left hand side. The max- 
imum value of k is of the order of the inverse of the lattice spacing, and 
in the thermodynamic limit the N vectors k tend to fill uniformly a d- 
dimensional sphere with a radius equal to this maximum value. Therefore, 
in this limit the possible values of k are distributed according to k d ^ 1 . It 
follows that in the thermodynamic limit only a vanishing fraction of the 
eigenvalues diverges (and at most as N, like Ao), also in the less favourable 
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case (concerning the distribution of the values of k near 0), when d — 1. 
This will be important for our computation of the partition function. It 
is also evident that negative eigenvalues are possible only for finite values 
of k, and that the least eigenvalue, in the case b = in (|lj]), is of order 
one in modulus. Thus, the all spectrum can be made positive by properly 
choosing the value of b in ([l9|). 

We can now estimate the behaviour of the rescaled eigenvalues per- 
taining to the interaction (||) . The eigenvalues of the rescaled interaction 
matrix are related to the Ak, being given by X^/N. If we choose a value 
of b such that the least eigenvalue Ak has a positive value e, then for the 
eigenvalues of the rescaled interaction the following relation holds: 

< 4 < ^ < 1. (24) 
N ~ N ~ 

According to what has been noted above, only a vanishing fraction will 
remain finite in the thermodynamic limit; this will be an important fact in 
the following. 

To end this subsection, we note that when a = the eigenvalues can 
be calculated from (p0|): 

^ = b + J2 e ~ ll " r = b-l+NS k0 . (25) 

r^O 

There are N — 1 eigenvalues equal to b — 1 and one eigenvalue equal to 
b + N — 1 = N. All the rescaled eigenvalues vanish in the thermodynamic 
limit except the single biggest one equal to 1. This is just the extreme case 
of the previous analysis. 



4.3 Canonical solution 

We now show that thermodynamic universality holds among all the rescaled 
a < d cases, in the sense that the specific free energy and the equation of 
state are the same. 

With the help of the analysis of the previous subsection, we are now able 
to find an exact solution of the rescaled model defined by the hamiltonian: 



1 N 1 S S N 



2N fr^ r\ 



with a < d, Jij — l/rfp N defined as in @2p, and = b defined so as to 
have (pi|). As in the mean field case, we have added an external magnetic 
field h. 
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In the following, whenever a site-dependent n-dimensional vector like Si 
occurs we indicate its components with double index quantities with the 
roman index varying in the range 1 ... N and the greek index in the range 
1 . . . n. The n-dimensional vector at site i is indicated with boldface letters 
Si while we let indicate the iV-dimensional row vector («Si M , . . . , Sn^), 
and S_ the corresponding transposed column vector. The components of 
a site- independent n-dimensional vector, like h in (|26|), will be indicated 
with, e.g., h^. 

The partition function for the model is given by eq. (||) with the new 
hamiltonian. Analogously to the mean field case, we define = fih^ and 
C — exp [—{P/2N)Y^ij Jij] = ex p(~-^/3/2); due to the site independence 
of the magnetic field, we have = [B^, . . . , B^). Introducing the matrix 
Rij = ((3/2N)Jij we can rewrite the partition function in matrix form: 

Z = cfd N @e^^ R ^ + ^\ (27) 

Like in the mean field case we want to make use of gaussian transformations 
to linearize the quadratic term. Then we first diagonalize the matrix R 
with the unitary matrix V such that VRV T = D, with Dy = RtSij, 
where Ri, the eigenvalues of R, are related to the eigenvalues Xi of Jij by 
Ri — (f]/2N)Xi. So we can write the first part of the exponent as 

£££3^ = (28) 

li in 

where = VS_„ . Because of ( pi| ) all the eigenvalues Ri are positive and we 
can apply the gaussian transformation (|J) to each term in the right-hand 
side of (pq), obtaining: 



[(4 7 r) Ar det R} 4 



HdviM-^^l , (29) 



where the appearance of det R in the denominator is due to the relation 
Yli Ri — det-R. Noting that R^ 1 = (D~ 1 )u, with the change of variables 
defined by = Vz^, and introducing the notation ]J ifi dzi^ = d Nn z, the 
previous expression can be written as: 

[(4tt) n dct R}9 J V ' 

Inserting this in @, and noting that £ p +^) = E l S i ' ( z * + B )> 

the partition function becomes: 



C 



[(Att) n det 



J d Nn ze -iy^R-% /■^q^s.k+b) (31) 
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This expression is similar to (||), but with an n-component integration 
variable z for each site i. With the same notations used in Section || and 
discussed in the Appendices, we therefore obtain: 



N 



d Ar„ 2e [-|E M ^ fl " 1 ^+E i lnG «-^l^+ B D] . (32) 



[(471-)" det J?]* 

This integral, analogous to (|o|) of the mean field case, can be solved with 
the saddle point method. However, here the justification requires more at- 
tention than before, since together with TV also the number of integration 
variables becomes very large. We postpone this point to the next subsec- 
tion, and for the moment we look for the stationary points of the exponent 
of ( Bp ) and we perform the stability analysis. We do not put in evidence a 
factor TV, since this is not necessary for the search of the stationary points 
and for the determination of their character. If we call /({zi^}) the expo- 
nent, the stationary points are given by the solutions of the system of Nn 
equations: 

fi = l,...,n i = l,...,N. (33) 



The value of the integral in (|32j ) will be determined by the absolute maxi- 
mum of f{{zifi})\ therefore, in the thermodynamic limit we will obtain for 
the specific free energy: 

1 

- 0F = lim — \nZ= -£ + lnfi„_i -nln2 (34) 

N— >oo IX 2 



max f In dct B In det [ — Hn 

2 2 V 2 



where Hq is the hessian of / computed in the absolute maximum (see below 
how it has to be interpreted in case of degeneracy of this maximum) . The 
essential difference with respect to the expression of the partition function 
of the a = case, Eq. (|i~2|), is represented by the terms containing det-ffo 
and deti?. There it was not necessary to consider the hessian of the ex- 
ponent; the reason is that ( |To| ) is n-dimensional, and the contribution of 
the hessian to the specific free energy vanishes in the thermodynamic limit. 
The second derivatives of f({z ifl }), necessary to compute the hessian, are 
given by: 

d 2 f 
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g ;_ 2 (K + Bl) ( ^ + ,^ )( ^+^ ) + (35) 



■ 9 - 2(|Zl + B|) ^TbF 

The system (|33|) is equivalent to the following: 



|z,+B|2 
+ g n - 2 (\z t +B\ 



N 



z,+B| 



Zip = 2y^J^jg n _ 2 (| 

/i=l,...,n z = l 



gj> + -Bp 
|z J+ B| 

. . . , N. 



B| 



(36) 



4.3.1 Homogeneous solutions 

Anticipating that the relevant stationary point is homogeneous on the lat- 
tice, we show what is obtained if we look for a solution in which z^ does 



not depend on i. Then (36) reduces to the system: 



^=/3 5n _ 2 (|z + B|)-^ 



B„ 



B 



,,n, 



(37) 



which is identical to (111): in fact, from the definition of the matrix Ra and 



from (|22|) we have that R-ij = f° r eacn i- We take the solution of 
( p7| ) already considered for a = 0, and we show that the stability analysis, 
now applied to the matrix (|35|), gives the same results as for that case. 
Afterwards, we will consider other possible solutions of ([56|), non homoge- 
neous on the lattice, and we will show that other possible maxima, among 
them, are only local. 

In the case with magnetic field we considered, without loss of generality, 
the vector h in the positive direction of the first axis, and we found that 
the corresponding solutions of ( ^7|) have = for /x ^ 1. For such zs, 
(and for = BS^i), the second derivatives vanish if ji ^ v, while for the 
others we have: 



d 2 f 



dzndzji 
d 2 f 



\(Rr 1 ) ij + 5 ij g , n _ 2 Qz 1 +B\) 



(38) 



9n-2(\Zl+B\) 

' \zi+B\ 



Therefore the matrix of second derivatives separates in n N x N blocks. 
For fi — 1 and for ^ 1 the eigenvalues of the corresponding matrix are 
given in terms of the eigenvalues Ri of Rij by, respectively: 



-~#r 1 +sn- 2 (i*i+£| 



1 

-2* 



g n _ 2 (\ Zl +B\) 
\zi + B\ 



(39) 
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From the definition of Rij and from ( pi| ) we have the following inequalities: 



Therefore, the same analysis presented in Appendix C, with equations 
and (p9|), can be performed, with the conclusion that the relevant stationary 
point is that with z parallel to h, while the other two solutions with z\ < 0, 
present for sufficiently small h and /3 > /3 C = n, are not maxima (except 
for n = 1, when one of the two is a local maximum). 

For h = 0, Eq. (|37|), as already seen, determines only the modulus 
z of the stationary point, becoming z — (3g n -i(z), which has always a 
solution z = 0, and, for (3 > (3 C , also a solution with positive z, which is 
infinitely degenerate for n > 1 and doubly degenerate for n = 1. In analogy 
with the a = case, one should then perform in (|32|), before applying the 
saddle point method, an angular integration over a global rotation. This 
would give in (|34| ) a further factor ln f?" , which does not contribute in the 



JV 

thermodynamic limit. In the point z — the hessian matrix is given by: 



d 2 f 



with eigenvalues: 



l^_r, , . „ „ 1 



= -~(^ + (41) 

z=0 2 n 

Using (pp[), we see that for (3 < (3 C this stationary point is a maximum, 
while for (3 > (3 C it is not a maximum. In the case of a point z* whose 
positive modulus z* satisfies the self-consistency equation for (3 > (3 C , we 
have: 



d 2 f 



-(R %• + Sij- 



- ^^[u,, A: l-^^V (.431 



dz^^dzjij 



We can look at this matrix as being given by the sum of different contri- 
butions, each one negative semi-definite. For the term multiplying <5 M „ this 
can be seen from (^0|), while for the term multiplying Sij we can refer to 
the analysis performed in Appendix C for expression (|96|). The degeneracy 
of the stationary point is reflected in the existence of n — 1 eigenvalues 
equal to zero; however, we have just noted that in this case an integra- 
tion over the region of degeneracy is implied before the application of the 
saddle point. This integration, that gives a vanishing contribution in the 
thermodynamic limit, is performed exactly over the directions correspond- 
ing to the eigenvalues equal to zero; as a consequence, these eigenvalues of 
det (— iJTo) must not be taken into account in the stability analysis. 
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4.3.2 Possible non homogeneous solutions 

Now we consider other possible solutions of (|3^) , non homogeneous on the 
lattice. We do not prove if they exist and if, in that case, they are maxima; 
rather, we prove that, if they exist and are maxima, they are local, i.e., the 
value of / is in those points smaller than for our homogeneous solution. 
Incidentally, we note that if a non homogeneous stationary point exists, 
then, because of translational invariance, all the configurations, obtained 
by the stationary one by any translation, are also stationary points with 
the same value. Let us begin by rewriting the stationary point equations 
( p6| ) in another form, posing z ifi + = w ifi : 

N 

w ifl = 2 J y~R ij g n ^ 2 (\^j\) T M 7+B^ 
~i l w il 

H=l,...,n i = l,...,N. (44) 

From this we can derive an equality, useful later, that is verified for our 
homogeneous solution. Taking, as before and without loss of generality 
= BSfji, and the relevant homogeneous solution with Wi^ = wS^i (i.e., 
Wi^ — for fi = 1 and wn — w > 0), we obtain: 

w = 0g n - 2 (w) + B, (45) 

where again we have used 2 R i j — f3 for each i . This is valid also for 
B = 0. Let us now consider another possible solution of ( f44| ) for which not 
all the |wi|s are equal. In this case we take the first axis in the direction 
of the v/i with the largest modulus, be it w/, and from ( f44| ) we have: 

N 

Wji 



w n = \w l \=2j2Ry9n-2(\^\)rA+Bi 

3=1 |W?I 

N 



3 = 1 
N 

< 2^i? y5ll _ 2 (|w / |) + S = /5 5 „_ 2 (|w i |) + B, (46) 

where we have used the monotonicity of g n -2{x) and that B\ < B. There- 
fore, a non homogeneous solution with different moduli for the w^s is such 
that all these moduli are smaller than that of the homogeneous solution 
satisfying (|45|). In fact, the properties of the functions g (see also Fig. 
1) imply that |w/| in ( [46| ) is smaller than w satisfying (ff5|). The same is 
true for a solution with all equal moduli but different directions; in this 
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case, in ( p6[ ) the first inequality becomes strict and the second becomes an 
equality. It is now sufficient , analogously to what we have done at the end 
of Appendix C for the a = case, to see that the exponent / in ( |32| ) is 
an increasing function of the moduli |w»|, and this will prove that all non 
homogeneous solutions of ( |36| ) are at most local maxima. It is not difficult 
to show, using (|44|), that, as a function of the w,-s, / in the stationary 
points is given by: 



^3n-2(|w l |)|w l | + ^5n- 2 (|w 4 |)^-^ + hi G n -2 (K |) 



(47) 



Posing Wj = XiSi, with X\ = |w{| and Si unitary vectors, ( J47| ) becomes: 

JV 



(48) 



The differentiation with respect to Xi gives, as it is easily seen: 

dj_ = 1^ _ g > n 2{xi ^ + lg' n _ 2 ( Xi ) Si . B. (49) 

The first term is always positive, according to property 5 of Appendix B, 
and also g' n -2( x i) i s always positive. This proves that for the homogeneous 
solution, satisfying ( fi"5| ) and for which • B = B for each i, f has the 
largest value. 



4.3.3 The free energy 

Now we are in the position to compute the terms in square brackets in 
( [34| ) and to write a more explicit expression for the specific free energy F 
in (p4j). We first write down an expression for max/. The maximum z* 
is homogeneous, then we have to know J^iji^ 1 )^- Since J2j ^hj = f 



for each i, then 



for each i. In fact, the first expression 



tells that an homogeneous vector is an eigenvector of R with eigenvalue ^ ; 
therefore, the same vector is an eigenvector of R^ 1 with eigenvalue %, that 
gives the second expression. It follows immediately that: 



max f 



N 



,*2 



2(3 



TnG„_ 2 (|z*+B|) 



(50) 



We now compute lndet (— \Ho)- Performing the stability analysis we 
have explicitly considered the eigenvalues of the hessian matrix in the homo- 
geneous maximum. Summarizing the results, we can write the expression 
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below for the eigenvalues of — \Hq, valid both when the maximum is at a 
positive value of z* (i.e., when B > and when B = but (3 > (3 C — n) and 
when it is at z* = 0. It is easily seen (always considering, for convenience, 
B along the first axis, and z* along the first axis also when B = and 
f3 > (3 C ) that the eigenvalues are given by: 

Ir^-P^z*) 

—R^ 1 —p2(z*) n—1 times (51) 

for i = 1, . . . , N. Here, when z* > 0, Pi(z*) = \g' n ^ z * + ^0 an d Vi{ z *) = 
Ja z*+b > wnue i when z* = 0, Pi(z*) = P2( z *) = ^- Therefore we have: 

N T /i \ / 1 \ n—1 



det = f[ (j^ 1 -PiCO) (J^ 1 -Pa( 

■ v Nn N 

- J (det i?)-" [] [(1 - 4i? lPl (z*)) (1 - 4i2 iPa (2*)) n - 1 ] . (52) 

' i=l 

The fact that we have to disregard the zero eigenvalues when B = and 
/3 > P c is automatically taken into account in the above expression. We 
then obtain: 

-i In det f-^#o^ = iVnln2 + ^fndet.R + 

N 

+ 5^pn(l-4J2 i p 1 («*)) + (n-l)ln(l-4R i p 2 («*))]. (53) 

i=l 

We will show in a moment that the last sum does not contribute to ([34j); 
then, using (|5(]) and (^3|) in (34), we have the following expression for the 
free energy: 

-0F = -| + lnn n _! - ^z* 2 + lnG„_ 2 (|z* + B|) (54) 

identical to ([12). We then obtain the same magnetization as for a — 
(see also fllOlJ) ), and the same equation of state (see ( fL4|) and (102)). 
Therefore, we still have a second order phase transition at /3 = j3 c , where 
M becomes 0. Of course we can make the same comments as for the mean 
field case (see the last paragraph of Section ||) , concerning the addition of 
the kinetic part. 

This concludes, apart from the remaining technical points treated be- 
low, our proof of the universality of all the models, for each n, when a < d. 



17 



We repeat that the difference, when n varies, is in the value of the critical 
temperature T c = — , but the overall behaviour is the same for all values of 
n. 

We are left with two points: the proof that the sum in (|5^) does not 
contribute to (Q), and the justification of the saddle point method. The 
next subsection is dedicated to the second point; here we treat the first. 
From ( f24| ) we have: 

(55) 

with the important specification that only a vanishing fraction of these 
eigenvalues remains finite in the thermodynamic limit. We can therefore 
write: 

1 N 

lim — V [In (1 - AR lPl (z*)) + (n - 1) In (1 - AR lP2 (z*))] = 

TV— *-OC iV L — ' 

i=l 

= ^ oo {-^E't 4 ^bi(^) + ("- 1 )^(^))] (56) 

+ ^ E " t ln (! - 4i W 2 *)) + (« - 1) ln C 1 - 4i2iPa(«*))]} , 

where the first sum is on the vanishing eigenvalues, and the second on the 
others; in the first we have substituted the logarithm with its first order 
approximation. We also point out that, following our stability analysis, the 
arguments of the logarithms in (|5^) are between and 1. We can indicate 
with rj the fraction of eigenvalues present in the second sum; if we also 
denote with R' max and R'^ax t ne largest eigenvalues in the two sums, the 
above positive expression is bounded from above by: 

(l-v)R' m ax + ri\H^-R'L x )\- (57) 

In the thermodynamic limit this expression goes to zero, since both R' max 
and r) go to zero. Therefore, the sum in ( p3[ ) does not contribute to the 
specific free energy (p3). 



4.4 Justification of the saddle point 

In the integral (^2|) the number of integration variables itself grows with 
N. The value of the maximum of the exponent diverges with N (this is 
seen in (jSCj)); but the use of the saddle point gives a proper evaluation of 
the integral if also the curvatures in all directions diverge with N, i.e., if 
all the eigenvalues of —^Hq diverge with N. We now show that this is not 
exactly verified, but since what we are interested in is the evaluation of the 
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specific free energy F (see (|34|) and (^4|)), we also show that nevertheless 
the saddle point method is still justified. 

The eigenvalues of — \Hq are given by (|5l]). The functions Pi{z*) and 
p%{z*) are finite, and for jR" 1 we can look at expression (HQ). According 



to our previous analysis in subsection 4.2, among the values of jR^ 1 only 
a fraction, vanishingly small in the thermodynamic limit, does not diverge 
with N. Along the directions corresponding to those eigenvalues the in- 
tegral should be computed explicitly, reserving the saddle point method 
to the other directions; however, we show that the error introduced using 
altogether the saddle point vanishes in the computation of F. 

Let us call collectively Si the eigenvalues of — ^Hq, where I runs from 1 
to P = nN , and indicate with v\, . . . , Up the integration variables in (|32|); 
then, the fact that the integral in that expression is evaluated with the 
saddle point is equivalent the the following replacement: 

d p ve f — ► e max/ / d p ve~^i=i S,v ' . (58) 



For the w/s for which Si does not diverge with N this replacement is not a 
good approximation, and we should more correctly write something like: 

e max/ /" dP^-EL^I-xK!-^) (59) 



with u > always, except that u = when vp +1 = . . . = vp = 0, and 
where 6p +1 , . . . ,Sp (with P < P) are the <5;s that do not diverge with N. 
The previous expression can also be written as: 



J dv P+1 . . .dv P e Si=p+i ^"i 

We are interested in the limit, when N — > oo, of -k times the logarithm of 
this expression. We have already seen that, when N — * oo, p ~ p goes to 
zero. This implies that the contribution of the last fraction in d60| ) to the 
evaluation of the free energy F vanishes in the thermodynamic limit, and 
this is equivalent to using the saddle point expression (pq). 



4.5 Microcanonical solution 

We briefly treat the point of the microcanonical solution of our system. The 
equivalence of different ensembles is a problem of general character, and 
here we only show that for our system the canonical and microcanonical 
ensembles are equivalent. 



19 



For a generic system, indicating collectively with T the coordinates of 
its phase space, its canonical partition function can be written as: 



Z{(3) = / dTe-^ = / dEu{E)e-P E = / dEe ^ B+u ^ , (61) 
J Jo Jo 



where uj(E) is the microcanonical density of states: 



u>(E) = J dTS(H(T) - E) . 



(62) 



In ( |6l| ) we have supposed that the hamiltonian is bounded from below at 
(as in our case) without loss of generality; the dependence on the number 
N of particles is not explicitly written. Introducing the specific energy 
U = § and using the definition of the microcanonical specific entropy in 
the thermodynamic limit 



5(f) = lim -]nu(E) 

TV— >oc iv 



(63) 



then, from the last expression in ( j6l| ) we see that the canonical partition 
function can be computed, in the thermodynamic limit, by the saddle point 
method, and the specific free energy is therefore given by: 

- 0F(P) = lim 1 In Z{f3) = max [-0U + 5(17)] . (64) 

N—>oo 1\ U 

If S(U) is convex, i.e., if 4tj| < 0, this relation defines a single value of U 



for each (3, U mc ((3), given by j§\ u=u 



— /3, and we can write: 
0F{0) = -f3U mc (p) + S{U mc {(3)). 



(65) 



On the other hand, if S(U) has a concavity region, as in the presence of 
first order phase transitions, then it can be easily deduced from ( |64| ) that 
the temperature derivative of F has a discontinuity. In our system we have 
only a second order transition, and in fact we have no discontinuity in ^ ; 
therefore 5(f) is convex and ( |65| ) holds (of course this is true whether or 
not we consider the trivial kinetic energy contribution to F). It is now easy 
to show the equivalence of the two ensembles. In fact the entropy 5 c (/3) 
computed from the canonical ensemble is obtained from (|65|) as: 



dF 



,8F 







— = f— = /r— 



8T 



dft 



Oft 



Umc(fi) - ^5(f mc (/3)) 



= (3 2 



da 



dp ^ SiUmcm 



1 dS_ 

~f3 dU 



dU mt 
u=u mc d P 



S(U mc (P)) 



(66) 
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using that j^j 



du I u=u, 



(3. Also for the canonical energy U{0) we obtain: 



U{(3) 




— {f}F{p)) = ^ 



S(U mc (f3))} 



dU rn 
~d(3 



m c 



U mc {(3). (67) 



5 Discussion and conclusions 

In this paper we have considered the problem of computing the partition 
function of lattice magnetic models with long-range couplings. We have 
studied a class of models in which the decay of the interaction with distance 
is gauged by the exponent a, smaller than the spatial dimension d in which 
the lattice is embedded. From the technical point of view, our study has 
required several steps: i) the analysis of the spectrum of the matrix R 
for a consistent application of the well known gaussian identity sometimes 
called Hubbard-Stratonovich transform (see (||) and (|30|)); ii) the analysis 
of a class of functions related to the modified Bcssel functions of the first 
kind, and characterized by the index n, the dimension of the spins; in) 
the application of the saddle point method to an integral with a diverging 
number of integration variables, which had to be justified; iv) the proof that 
possible stationary points, if any, in the general case (in the strict sense 
< a < d) , of the exponent in ( |32| ) can be at most only local maxima, and 
are therefore irrelevant in the thermodynamic limit. 

In our computations, we have not explicitly considered the case a — d 
and we have restricted the long-range couplings to a power form. However, 
it is not difficult to argue that also for a power decay with a = d and 
for a more generic form of Jij with a long-range character (i.e., X^Li \ 
diverging with N) we would have obtained the same results. In fact, the 
basic points for our computations are: firstly, the divergence with N of 
the quantity S in ( |2l| ) and consequently of N in (pffi), and secondly the 
fact that only a vanishing fraction, in the thermodynamic limit, of the 
eigenvalues Ak diverges, with the consequence, given the first point, that 
only a vaninshing fraction of (see (§4J)) remains finite. When a = d 
all these points would remain; the behaviour of the quantity S in (|2l]) for 
large N would then be given by IniV. Also for a generic Jy all these points 
would still be true, even if the divergence law of S with could possibly 
be difficult to write explicitly. 

Let us remark that we have also shown the equivalence, for these models, 
of the microcanonical and the canonical ensembles, in spite of the lack of 
additivity. 

In conclusion, we have found an entire class of lattice spin models with 
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a universal thermodynamic behaviour. 

These results on lattice systems can be helpful for the thermodynamics 
of continuous long-range systems; it should be stressed, however, that, 
without the Kac's prescription, not only additivity, but also extensivity is 
violated, and this could imply ensemble inequivalence. The microcanonical 
ensemble would then be the natural framework for the study of those cases. 

Appendix A: The functions G p (x) and g p (x) 

In (^) we have, for each site i, an integral of the form: 

s 



d& e & z (68) 

extended on the surface of the unit sphere in dimension n > 2, where the 
unit vector S lies. For n = 1 the integral is substituted by the sum over 
S = ±1; in that case also z has a single component z, ( |68| ) is given by 
2cosh(z), and we do not need further analysis. For n > 2 we can choose 
the axes of the unit sphere such that z lies along one of them; besides, we 
can introduce polar coordinates in n dimensions, taking as the polar axis 
the one along z. Then, if z = |z|, it is easy to show that ( fig) ) becomes: 



d0 e 8- * = n n _i / desm n - 2 6e zcose , (69) 



where 9 is the polar angle, while Q n is the area of the unit sphere in n 
dimensions (with f2i = 2); it can be expressed in terms of the gamma 
function as O n = 27r"/ 2 /r(ri/2). This last factor is not considered further, 
and for the integral in (|6^) we introduce the following notation: 

G p (x)= d6sm p 9e xcosf> = 2 / d9 sin p 6 cosh(x cos 6), (70) 
Jo Jo 

where the parameter p is related to the dimension n > 2 of the spin vector 
by p = n — 2, and therefore p takes the non neg ative integer values. If 
n = 1 it is understood in Section |^ and Section [f.j that G-\(x) = cosh(x) 
and f^o = 2. It can be shown (see next appendix) that G p {x) is related to 
the modified Bessel function of first kind of order |, i.e. 

Ir (x) 

G p (x) oc -ai-i, (71) 

£2 
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relation that is valid also for p = — 1. The other functions of interest are 
the derivative of G p (x): 

G' p {x) = d9 cos 8 sin p 8e x cos 6 = 2 / d8 cos 6 sin p sinh(x cos 8) 
Jo Jo 

(72) 

and the logarithmic derivative: 

, . d . . G;(x) C^cosflsin^e^ 08 * 

Jp(i) = t lnGp(x) - - 



pw G p (ar) /* d0 sin p 6e x cos 6 

JJ r/ ' 2 d9 cos 6* sin p sinh(x cos ( 



J ( ^ 2 d8 s'm p 6 cosh(x cos 8) 



(73) 



For n — 1 (i.e., p = — 1) they are substituted by sinh(x) and by tanh(x), 
respectively. In the next appendix we prove some properties of the functions 
G p (x), G' p (x) and g v {x), which are necessary for the analysis of the self- 
consistency equations ([TJ]) or (|37|). 

Appendix B: The properties of g p (x) 



As we have seen, in (|70|), ( |T2[ ) and (|7^) the parameter p takes the non 
negative integer values. We note however that, if we take p to be any real 
number, the functions inside the integrals are integrable for any p greater 
than —1, and besides it can be easily proved that the limit for p — > —1 of 
( |73|) is exactly tanh(a^). As a consequence all the properties that we will 
show for G p (x), G' p (x) and g p (x) are valid for any p > —1 (if necessary, for 
p = — 1 the explicit expressions that we have given can be invoked). Some 
properties will be immediately evident from the definitions, while others, 
especially those of g p (x), will need a proof. 

Considering first G p (x) and G' p (x), we begin with their properties of 
symmetry with respect to x inversion and monotonicity. In particular, 
the positive function G p (x) is even in x and is monotonically increasing 
for x > 0. This is most easily seen from the second expression in (f70|). 
From the second expression in J72] ) it is also readily seen that the odd 
function G' p (x), positive for x > and negative for x < 0, is monotonically 
increasing. A relation which is useful for the following can be derived with 
an integration by parts in fl72j), namely: 



G' p (x) = —-G^x). (74) 

We also give the value of the first few derivatives of G p (x) in x = 0. 
Since G p {x) is even, the odd derivatives vanish for x = 0, while from the 
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expression 

G p m \x) = -— G p {x) = / d6 cos rn 8 sm p 0e xcos 9 (75) 
dx m J 

integrations by parts, after posing x = 0, give, for the second and fourth 
derivatives in x — 0, the expressions: 

G?(0)=^Gp(0), ^ ) (0)= (p + 4 + 2) ^(0). (76) 

These are useful for the computations of the necessary derivatives of g p {x) 
in a; = 0). Before turning to g p (x), it is interesting to derive a second 
order differential equation satisfied by G p (x). From the relation G p (x) = 
GJx) — G p+2 {x), obtainable from ( f75| ) posing cos 2 8 — 1 — sin 2 0, using 
(l?4j) we have: 

x 2 G' p \x) + (p+ l)xG' p {x) - x 2 G p (x) = 0. (77) 
From this it is easy to derive the equation satisfied by W p (x) = x*G p (x): 



x 2 W' p \x) + xW p {x) 



x 2 



(I)" 



W p (x) = 0, (78) 



which is the modified Bessel equation with parameter | . Knowing also the 
limiting properties, for x — * 0, of G p (x), it can be concluded that W p (x) is 
proportional to the modified Bessel function of first kind of order | , Ip (x) . 
Regarding g p {x), we show the following properties: 

1. g P (x) is odd, positive for x > and negative for x < 0; 

2. g p (x) is, in modulus, smaller than 1, and it tends to ±1 for x — ► ±oo; 

3. g' p (x), which is even, is positive; 

4. 9p ^ ; which is even and positive, has a derivative which is negative 
for x > and positive for x < 0; 

5. SpIeI — g' {x), which is even, is positive; 

6. g'p(x), which is odd, is negative for x > and positive for x < 0. 

We give the proofs of these properties. 

Property 1. It is an immediate consequence of the definition and of the 
simmetry and monotonicity properties of G p (x) and G' p (x). 

Although also the second can be easily derived from the definition, 
we prove it, together with the remaining, by introducing the first order 
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differential equation satisfied by g p {x). From the definition of g p {x) we 
have: 



G p (x) \Gp(x) J Gp(x) p 
Substituting G p (x) as a function of G' p (x) and G p (x) from (|77[), we have: 

g> p ( X ) = l-P±±g p ( X )-gl( X ). (80) 

This equation must be supplied with the initial condition g_p(0) — 0. The 
first few derivatives of g P (x) in x = can be obtained from ( p0| ) or from the 
derivatives of G p {x). In particular, since g p (x) is odd, the even derivatives 
in x = vanish, while for the first and third derivatives we have: 

4<°> = <' (0) = - (p + 4)( P + 2)2 - ( 81 ) 

We now proceed with the remaining proofs, using (|o|) and ([Si]). 

Property 2, first part. g p (x) starts at x = with positive derivative 
(see ([si])). Therefore, if, for x > 0, it intersects the value g p = 1, it must 
do that at least once with non negative derivative. But from (^0|) we get 
for g p — 1 that g' p — — ^r-j which is negative for finite x > 0. Then, g p (x) 
remains, for x > 0, between and 1. 

Property 3. By differentiating (|8(]), we get the following equation for 

9» = P -^T9p{x) P -^9' P {x) - 2g p (x)g' p (x). (82) 

Since ^(0) = ^5 > ^ ; if g' p {x) intersects the value 0, it must do that at 
least once with g' p ' non positive for x > and non negative for x < 0. But 
from ( p2| ) we get for = that g' p ' = ^r-g p , which is positive for finite 
x > and negative for finite x < 0. Then, g' p (x) > for all x. 

Property 2, second part. We restrict to positive x; the proof for negative 
x follows from the fact that g p (x) is odd. Since g' p (x) > and g P (x) < 1 for 
finite x, g p (x) must have an asymptote, and for x — > 00 g' p (x) must tend 
to 0. In that limit ( |S0| ) becomes = 1—^. Therefore, for a; — > 00 g P (x) 
tends to 1. 

Property 4. We define h p {x) = From @ it is easy to derive the 

differential equation satisfied by h p {x): 

h'Jx) = -- P —h p {x) - xh 2 p (x). (83) 

^ X X 1 

Following what has already been proved, h p (x) is even and positive, and 
we have: 

M0) = ^, h>(0)=0, K(0) = - (p + A) l + 2r (84) 
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By differentiating (|3|) we have the equation: 

h£(x) = ~^ + ^T h p( x ) - ^rKW ~ h l( x ) ~ 2xh p (x)h' p (x). (85) 

From (|84| ) we see that, for x positive sufficiently small, h' p (x) is negative; 
therefore, if for larger x h' p (x) becomes equal to 0, it must do that with a 
non negative h' p \x), and with h p (x) smaller than — But from (^5|) we 

have for h'Jx) = that h p (x) — —h p (x) — [1 — (p + 2)h p (x)], which is 
negative. Then, for x > h' p (x) is always negative, and it can also be easily 
seen that it tends to zero, together with h p (x), for x — > oo. The proof for 
negative x follows from the fact that h p {x) is even. 

Property 5. This is a simple consequence of the previous property, once 

we note that hp {x) = 9 -^- 9 -^ = ± (g' p { x ) - S&l) < 0. 

Property 6. By differentiating (|8^), we get the following equation for 

<(*) = -2^ 5p (x) + 2^±i. 9 ;^) + 

^±l ff ;'(^)-2 ff ; 2 (x)- 2 ffp (a;)^(a;). (86) 

The positive function <7 p (cc) starts at x = with the value and for suf- 
ficiently small x g p (x) is negative (see ( 
becomes equal to 0, it must do that with a non negative g p '(x). But from 

© we have for g%(x) = that 9p "(x) = -2g' p \x) - 2^+1 (^i _ g ^ , 
which is negative, using property 5. Thus, for x > 0, g p {x) is negative. 
The proof for negative x follows from the fact that g p (x) is odd. 



81)); therefore, if for larger x g p (x) 



Appendix C: Mean field solution of the a = 

case 

We give here some details of the mean field solution. We start from the 
stationary point equations (pi]), that we rewrite here: 

f)f 1 7 + P, 

-L = - ^+g n _ 2 (\ z + B\)^-£ = } M =l,2,...,n. (87) 

Let us begin with the case h > (we remind that B — (3h). Without loss 
of generality, we can take the magnetic field h in the positive direction of 
the first axis, i.e., h\ = h and = for fx ^ 1. Then, from (|87| ) we have, 
for [i 7^ 1: 

z M = /?3n-2(|z + B|)— J— , M = 2,...,n. (88) 
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The possible solutions of these n — 1 equations are: z^ = for each jj, = 
2, . . . , n, or, if some is not zero, such that /3 g "^j^ B ^ = 1. But in the 
second case the equation for fi = 1 would become z\ = Z\ + B, which is 
not acceptable since B > 0. Thus, = for fi > 1, and the remaining 
equation for /i = 1 becomes: 

Zl =ftfa-a(l*i + #1) i^ 1 (89) 

From the symmetry properties of g n (x), described in the previous Ap- 
pendix, this equation is equivalent to: 

zi = 0g n -a(zi + B). (90) 

This equation can be solved graphically. We use the properties of the 
functions g n {x) that we have proved, and for a visual aid we refer to Fig. 
1. In particular, since the maximum of the positive function g' n _2{x) is for 
x = and is equal to ^, and the odd function g'n-2{ x ) * s negative for x > 0, 
we have that (|9^) always admits a single positive solution for Z\, while, if 
/3 > /3 C = n, it can also have, for sufficiently small h, two other negative 
solutions. From the stability analysis we will see that the relevant solution 
is the positive one. 

We now turn to the case h — 0. The stationary point equations become: 



dz~ 



Zfi + g n -2{z)— = 0, fi = 1,2, ...,n (91) 



with z = |z|. It is readily seen that these equation determine only the 
modulus z. In fact, posing z^ — zSu, we see that the unit vector s, giving 
the direction of z, is left free, and the equation for z is: 

z = Pg n - 2 {z), (92) 

which is of the same form of (|90|). Again, from the properties of g n {x) and 
looking at Fig. 1, we have that for j3 < j3 c = n the only solution is z = 0, 
while for (3 > (3 C we also have a positive solution. The stability analysis 
will show that in this last case the relevant solution is the positive one. 
To summarize, in both cases, with or without magnetic field, the relevant 
stationary point z* is such that its modulus z* satisfies the self-consistency 
equation (^), with the further characteristics specified soon after it. 

Let us now show the stability analysis. The second order derivatives of 
/(z) are given by: 

- -^ AW +g n _ 2 (|z + B|)— £ — ,^ D|2 (93) 



dz^dzu P » n -^ u |z + B| 
5 „_ 2 (|z + Bl) (^ + B,)(z +B^ ^ ^ lz + Bl) 



BP 1 "/••'-"--vi- 1 ""\ Z + B\ 
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They form the hessian matrix, that, at the stationary point, is negative 
definite if this point is a maximum. 

Let us begin with the case without magnetic field, h = 0. In this case 
the second derivatives become: 

^ / 1 r . / I \ z n^v 

dz»dz v - /3 V+S "- 2W z 2 

~ 9n-l{z)^- + Saug n -2{z)-. (94) 

z 6 z 

We have shown before that there is always a solution z = 0. Using the 
properties of the function g n (x), the hessian matrix at this point is given 

by 



d 2 f 
dz^dz^ 



z=0 

and is already in diagonal form. We immediately see that above the critical 
temperature (i.e., for (3 < (3 C = n) the only stationary point z = is a 
maximum; for (3 > (3 C it is a minimum and therefore it is not the right 
solution. At the critical temperature [3 — n the second derivatives vanish 
and the saddle point in (|l0|) must be solved considering derivatives of higher 
order, showing again that z = is a maximum. We do not show it here. 
For the solution with positive z for (3 > (3 C (that has an undetermined 
direction if n > 1 and that we call z* ) , using the stationary point equation 
"32|) we obtain for the hessian matrix: 

. 9 ; l - 2 (^)-^4 f!1 V 06) 



d 2 f 
dz^dzy 

The term in brackets is negative, according to the properties of the function 
g n {x), and we only have to study the matrix A M „ = z*z*. It is very easy to 
show that this matrix has an eigenvalue equal to z* 2 and n — 1 eigenvalues 
equal to 0. Therefore, coherently with the fact that for h = the stationary 
point is infinitely degenerate if n > 1, we find that this point is a maximum, 
but if n > 1 there are n — 1 directions along which /(z) does not change. 

We now consider the case with a magnetic field. We have studied the 
stationary points with hi — h and — for fi =/= 1, for which z^ = for 
fi 7^ 1. In general, for any z such that z M — for /j ^ 1 (and = £><5 M i), 
the hessian matrix becomes: 



d 2 f 



dz fl dz u 



(97) 
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which is again diagonal. If we consider first the diagonal term with ji = 1, 
we see that it is equal to 



a 2 / 

dz 



z,,=zi<5 fl i 



9 n -l(\z\ + B \ 



(98) 



The solution with z\ > is always present, and it is clear (see Fig. 1), 
that in that case (|9^) is negative. Among the two solutions with z\ < 0, 
which exist if > (3 C = n and if h is sufficiently small, it is also clear 
that for the one with the smaller \zi\ (|98| ) is positive, while for the one 
with the larger \zi\ ( p8| ) is negative (for the particular value of h for which 
these two solutions coincide (|98| ) is zero). Thus, if n = 1, also one solution 
with negative Z\ is a maximum; it is the metastable solution, with the 
magnetization opposite to the magnetic field, found in the hysteresis curve. 
We will see below its metastability, i.e., that for such Z\ the exponent in 
( |To|) is a local maximum, smaller than the absolute maximum obtained for 
the positive solution. Now we show that when n > 1 the only maximum is 
that with the positive z\. In fact, for /i^lwe have, for z\ satisfying (Bw: 



dzl 



B 



(99) 



These terms are negative for the solution with positive Z\, while they are 
positive for the solutions with negative z\ (for which z\ < z\ + B < 0). 

At the end of this computation we can say that in the thermodynamic 
limit the specific free energy F, obtainable from (|l~0|), is given by (|l2|), that 
is rewritten here: 



0F = 



lim — In Z 

N^oo N 







In fi n _i + max 



2(3 



z 2 +lnG„_ 2 (|z + B|) 



= -^ + lna i -i-^z* 2 + lnG„- 2 (|z*+B|). (100) 

2 2fj 



The magnetization is simply given, using (pl\), by 

8 



1 d 

M = lim ^—\nZ 

n^oo N/3 dh 

dz* 



E 
E 



i 



(3 



+ 3 „_ 2 (|z* +B| 



ff„_ 2 (|z* +B|) 



Iz*+Br 



^ + Bp 



B| V dB 



(101) 
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where — x M is the unit vector in the direction of the /i-th axis. As 
explained in the main text, while the degeneracy of z* when h = and 
(3 > (3 C does not affect ( |100D , the actual direction of M in a real system in 
this case is determined by the boundary conditions. The equation of state 
is obtained computing the specific energy, using (|7|) and ( |l()l| ): 

1 d 

U = — lim — — — In Z = 

n^oo N d/3 

1 z* 2 z* dz* , z*+B fdz* 

~ 2~W + J'W~ 9n ~ 2 (|z + l} FTb] ' \W 

= 1 --^-**. h= 1 - { l-M>)-hM, (102) 

2 2/3 2 (3 2 V ; v ' 

where M = |M|. 

We end this Appendix by showing that for n = 1 the maximum of the 
exponent in (^) with the magnetization opposite to the magnetic field is 
only local. We do this with some relations that are valid for any n, since 
they are useful for the study of the general model. Using (p7|), the exponent 
in ( p~0| ) in any stationary point, when this is in the same direction of the 
magnetic field (assumed to be the positive direction of the first axis), is 
given by: 

1 1 xB 

- -3n-a(M) • M + ^n-adsD-j-r + lnG„_ 3 (|s|), (103) 

2 2 \x\ 

where x = Z\ + B; it can also be written as: 

- -g n _ 2 {\x\) ■ \x\ + -g n - 2 (\x\)Bsign(x) + lnG n _ 2 (|x|). (104) 

The derivative of this expression with respect to |x| gives: 
1, , ( 9n-2(\x\) , .A , I , 



2 N { \x\ ~ g 'n-^) + ^ n -M)Bsign{x). (105) 

According to the properties of the function g n (x), the first term and the 
coefficient of sign(x) in the second term are positive. Knowing also that the 
solution with magnetization opposite to the magnetic field has a value of |x| 
smaller than that of the solution with h and M in the same direction, this 
is sufficient to prove that the former, for n = 1, is only a local maximum 
of the exponent in (|l0|). 
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Figure 1: In both panels the curved lines represent the qualitative be- 
haviour, for any n and for a range of z about 0, of f3g n {z + B) for /3 < (3 C 
(solid line) and for f3 > (3 C (dashed line). The straight line is the bisectrix, 
a) B — 0: when (3 > (3 C the relevant solution for the modulus z in (92) is 
the positive intersection, and not z = 0; b) B > 0: the relevant solution 
for zi in (|90|) is always the intersection with Z\ > 0, even when j3 > f3 c . 
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